Draft version February 2, 2008 

Preprint typeset using L^T^ style emulateapj v. 14/09/00 



DYNAMICAL EVOLUTION OF NEUTRINO-COOLED ACCRETION DISKS: DETAILED 
MICROPHYSICS, LEPTON-DRIVEN CONVECTION, AND GLOBAL ENERGETICS 

William H. Lee', Enrico Ramirez-Ruiz^-^ and Dany Page' 

Draft version February 2, 2008 

ABSTRACT 

We present a detailed, two dimensional numerical study of the microphysical conditions and dynamical evolu- 
tion of accretion disks around black holes when neutrino emission is the main source of cooling. Such structures 
are likely to form after the gravitational collapse of massive rotating stellar cores, or the coalescence of two com- 
pact objects in a binary (e.g., the Hulse-Taylor system). The physical composition is determined self consistently 
, by considering two regimes: neutrino-opaque and neutrino-transparent, with a detailed equation of state which 

' takes into account neutronization, nuclear statistical equilibrium of a gas of free nucleons and alpha particles, 

^SJ I blackbody radiation and a relativistic Fermi gas of arbitrary degeneracy. Various neutrino emission processes are 

considered, with capture onto free nucleons providing the dominant contribution to the cooling rate. We find 
^ ■ that important temporal and spatial scales, related to the optically thin/optically thick transition are present in the 

disk, and manifest themselves clearly in the energy output in neutrinos. This transition produces an inversion of 
the lepton gradient in the innermost regions of the flow which drives convective motions, and affects the density 
and disk scale height radial profiles. The electron fraction remains low in the region close to the black hole, and 
if preserved in an outflow, could give rise to heavy element nucleosynthesis. Our specific initial conditions arise 
from the binary merger context, and so we explore the implications of our results for the production of gamma ray 
bursts. 

■ Subject headings: accretion — disks — dense matter — gamma rays: bursts — hydrodynamics — neutrinos 
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1. INTRODUCTION 



O 

lO . Gas accretion by a concentration of mass is an efficient way to transform gravitational binding energy into radiation (Salpeter 
' 1964; Zel'dovich 1964), and is responsible for observable phenomena at all scales in astrophysics. The more compact the object, the 
greater the efficiency. In many cases, systems are observed in steady or quasi-steady state during this process (e.g, CVs, LMXBs 
O i' and AGNs). Some external agent (e.g., the interstellar medium, a companion star) provides a continuous supply of mass, energy 
and angular momentum over a timescale that is much longer than the accretion time. The energy dissipated by the flow before 
it is accreted, and more importantly, what happens to this energy, determines many of the properties of the flow itself. In most 
^ cases cooling is present through some radiative electromagnetic process, which acts as a sink, removing the dissipated energy. In 
, the standard thin-disk theory developed by Shakura & Sunyaev (1973), it is extremely efficient, maintaining a "cool" (in the sense 
that kT ^ GMirip/r, where M is the mass of the central object) and thin (with a scale height H ^ r), nearly Keplerian disk. The 
requirements of steady state, hydrostatic equilibrium in the vertical direction, energy and angular momentum balance, and an equation 
K> of state yield a solution for all relevant variables as a function of the disk radius, r. The key question of angular momentum transport 
5_j was addressed by Shakura & Sunyaev with their famous a prescription, which allows a parametrization of the viscous stresses and 
Ci ^ energy dissipation rates. Magnetic fields could be at the physical origin this effect, through the magneto-rotational instability (MRI, 
Balbus & Hawley 1991; Hawley & Balbus 1991). In the case of inefficient cooling, a solution was found in the 1970s for so-called 
slim disks (Shapiro, Lightman & Eardley 1976), and later for a class that became known as ADAFs (Advection Dominated Accretion 
Flows, Ichimaru 1977; Narayan & Yi 1994, 1995; Abramowicz et al. 1995). In this case, the cooling is negligible, either because 
the optical depth is large, or because the density is so low that the radiative efficiency is extremely small. The dissipated energy is 
advected with the flow and may be absorbed by the central object. Such flows are geometrically thick, with scale heights // ~ r. 
What little radiation does emerge from them arises, as in the thin disks, through various electromagnetic processes. 

In either of the above scenarios, accretion (and the accompanying radiation) is usually thought to be limited by the Eddington rate, 
a self-regulatory balance imposed by Newtonian gravity and radiation pressure. The standard argument gives a maximum luminosity 
^Edd = 1.3 X IQ^^(M /Mq) erg s"'. Although this may not be strictly the case in reality — as in the current argument concerning the 
nature of ULXs observationally (Rappaport, Podsiadlowski & Pfahl 2005; King & Dehnen 2005), and also quite general theoretical 
considerations (Abramowicz 2004) — it does exhibit the qualitative nature of the effect of radiation pressure on accreting plasma, in 
the limit of large optical depth. 

The problem is circumvented (or at least deferred by nearly sixteen orders of magnitude in luminosity) if the main cooling agent is 
emission of neutrinos, instead of photons. This regime requires correspondingly large accretion rates, of the order of one solar mass 
per second, and is termed hypercritical accretion (Chevalier 1989). It is important, for example, in the context of post-supernova 
fallback accretion onto a proto-neutron star. In such a situation, the densities and temperatures are so large {p ~ 10'^ g cm"^, 
T ~ lO" K) that photons are completely trapped, and energetic neutrinos are emitted in large amounts, cooling the gas and allowing 
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accretion to proceed. The reader may wish to consult the very clear introduction in this context given by Houck & Chevalier (1991). 
When the stellar envelope experiences fallback onto the central object, the system may be considered to be nearly in a steady state, 
because of the timescales involved. There are other, more violent situations, relevant to the study of gamma-ray bursts (GRBs) in 
which the assumption of steady state is not justified because of a time-varying mass and energy supply, and which therefore require 
a dynamical analysis for their proper description. 

In general, whether the system cools via electromagnetic radiation or neutrinos, analytic steady state solutions are found by 
assuming that the cooling is either efficient (thin disk) or inefficient (thick disk). Over the last several years, many groups have 
considered these cases in the neutrino cooled regime (Popham, Woosley & Fryer 1999; Narayan, Piran & Kumar 2001; Kohri & 
Mineshige 2002; DiMatteo, Perna & Narayan 2002; Yokosawa et al. 2004) and established the general features of the solutions. 
The real solution may be very different, in particular (and in what concerns us for the following) because of the previous history of 
the gas that constituted the accretion disk and how the disk formed. It turns out this may be quite important in the case of GRBs. 
Numerical analysis of the evolution of the disk with no assumptions concerning the steady state has been carried out by Setiawan et 
al. (2004). They have reported a three dimensional calculation with detailed microphysics. Unfortunately, computational limitations 
do not allow one to continue such calculations for more than approximately 50 ms, whereas the typical short ORB lasts 0.2 s. 

The sources of GRBs are now established (at least based on observations of X-ray, optical and radio counterparts) to lie at 
cosmological distances, with redshifts z ~ 1 - 4 (van Paradijs, Kouvehotou & Wijers 2000, and references therein). At such scales, 
the absolute energetics of each event is approximately 10^^ erg, assuming isotropic emission. Apparent coUimation inferred from 
achromatic breaks in the afterglow light curve may reduce this to 10^" erg, depending on the particular event (Frail et al. 2001; 
Panaitescu & Kumar 2001; Berger, Kulkarni & Frail 2003). Although one generally considers two classes of GRBs based on 
duration (shorter and longer than about 2 s, Kouveliotou et al. 1993), all afterglows to date (and the corresponding inferences) come 
from long bursts (although see Lazzati, Ramirez-Ruiz & GhiseUini 2001). Strong evidence in favor of a SN/GRB association 
(Kulkarni et al. 1998; Stanek et al. 2003; Hjorth et al. 2003) now comes from several events (albeit all relatively nearby), hinting 
at an underlying physical connection between the two phenomena. This is the central basis of the collapsar model (Woosley 1993; 
MacFadyen & Woosley 1999), in which rotation of the pre-supernova star allows for the creation of a massive accretion disk, fed by 
the infaUing stellar envelope. The GRB is then powered for a fallback time, which can be long enough to account for the observed 
durations. In the case of short bursts there is less evidence to go on, but one possibility is that they arise from compact binary mergers, 
with various combinations of black holes and neutron stars in the binary (Lattimer & Schramm 1976; Paczynski 1986; Eichler et al. 
1989; Paczynski 1991 ; Narayan, Paczynski & Piran 1992) such as in the first binary pulsar to be discovered (Hulse & Taylor 1975). 
The tidal disruption of one star by the other gives rise to an accretion disk which powers the GRB. A black hole is either present 
from the start (as a variation of the Hulse-Taylor system) or is produced by the collapse of a supramassive neutron star shortly after 
the merger itself. In this case there is no external agent feeding the accretion disk, and thus the event is over roughly on an accretion 
limescale (which would be on the order of one second). Clearly, investigating either of these scenarios requires time-dependent, 
multidimensional calculations with detailed microphysics (Ruffert, Janka & Schafer 1996; Kluzniak & Lee 1998; Ruffert & Janka 
1999; Lee 2001; Rosswog & Ramirez-Ruiz 2002; Rosswog, Ramirez-Ruiz & Davies 2003; Rosswog, Speith, & Wynn 2004). The 
energy released by accretion is then transformed into a relativistic outflow which produces the gamma ray burst (see Meszaros 2002; 
Zhang & Meszaros 2004; Piran 2004, for reviews). We note that either scenario would produce a distinct gravitational wave signal, 
which could in principle be detected by interferometric systems such as LIGO in the future, thus establishing the nature of the 
progenitor system without a doubt. 

The regime in which the gas lies is unUke any other commonly encountered in astrophysics. The high densities and temperatures 
lead to photodisintegration of the nuclei and the establishment of nuclear statistical equilibrium (NSE). Furthermore, neutronization 
becomes important in the innermost regions of the flow and weak interactions determine the composition, with the electron fraction 
falling substantially below 1/2, and the gas correspondingly becoming neutron-rich. If this composition is somehow frozen and 
transported out of the gas and into an outflow, interesting nucleosynthesis of heavy elements could occur (Qian & Woosley 1996; 
Pruet, Woosley & Hoffman 2003; Pruet, Thompson & Hoffman 2004). Realistic physics input of this kind allow us to obtain more 
reliable estimates of the actual energy released from the disk, and potentially available to power a GRB. An added complication, 
which affects the composition (Beloborodov 2003) is that even if the main cooling mechanism is neutrino emission, these are not 
entirely free to leave the system, as scattering (mainly off free nucleons) is important enough to suppress the emission in the dense 
irmer disk. We find that in fact the opaqueness of the material may lead to convection through the establishment of a composition 
gradient that does not satisfy the classical requirements for stabiUty. To our knowledge, this is the first time that this has been 
addressed in this context. 

In previous work we initially studied the merger process for black hole-neutron star binaries in three dimensions, paying particular 
attention to the structure of the accretion disks that formed as a result of the tidal disruption (Lee 2001). Follow-up work used the 
results of these simulations, mapped to two dimensions in azimuthal symmetry, as initial conditions with simple (ideal gas) input 
physics, and no realistic cooling included (Lee & Ramirez-Ruiz 2002). A first approximation at a reaUstic equation of state and the 
effects of neutrino opacities was reported more recently (Lee, Ramirez-Ruiz & Page 2004). 

In this paper we improve upon our earlier results in three important ways, that make them more realistic, and particularly relevant 
to the study of GRBs. First, we use a much more detailed equation of state, appropriate for the actual physical conditions found 
in post-merger accretion disks. This includes a relativistic electron gas of arbitrary degeneracy, radiation pressure and an ideal 
gas of free nucleons and a particles. The composition is determined self-consistently by considering weak interactions in two 
different regimes: neutrino-opaque and neutrino-transparent. Second, we include neutrino emission as the main source of coohng, 
by considering the relevant reaction rates (electron and positron capture onto free nucleons, bremsstrahlung, pair annihilation and 
plasmon decays), taken from tables and fitting formulae valid over wide ranges of temperature and density whenever possible. Third, 
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we compute an approximate optical depth for the fluid to neutrinos, using scattering off free nucleons and a particles as a source of 
opacity. The emission rates and pressure due to neutrinos are then suppressed and enhanced, respectively, by an appropriate factor. 

We begin this paper by a presentation of the physical conditions Ukely to occur in post-merger accretion disks in §2, immediately 
following the coalescence of the two compact objects. Our results follow in §3, and a discussion of these, and the implications they 
have for GRBs is presented in §4. 

2. POST-MERGER ACCRETION DISKS 

The accurate study of the dynamical evolution of the accretion disks requires knowledge of the physical conditions within them, 
and the use of an appropriate equation of state. Below we give a general overview of the physical conditions in the disk, followed by 
a detailed presentation of the equation of state and the effects of neutrinos. 

2. 1 . Physical conditions 

The accretion structures that are formed as a result of the merger of two neutron stars or the tidal disruption of a neutron star 
by the black hole become azimuthally symmetric fairly quickly, within a few tens of milliseconds (Ruffert, Janka & Schafer 1996; 
Lee 2001). These disks typically contain a few tenths of a solar mass, and are small, with the bulk of the mass being contained 
within 200 km of the black hole (which harbors about 3-5 solar masses). The disks are dense (lO'' < p[g cm~^] < 10'~), with high 
internal energies (10'° < r[K] < lO") (see e.g., Ruffert & Janka 1999; Rosswog et al. 1999; Lee & Ramirez-Ruiz 2002). In fact, 
the temperature is high enough that nuclei become photodisintegrated and there is a mixture of a particles, free neutrons and protons, 
electrons and positrons. The timescale for /3-equilibrium is given by f/j w (aNenj^c)~^ , where a^e is the cross section for capture 
onto free nucleons and n^ is the number density of free nucleons (see e.g., Shapiro & Teukolsky 1983). Under these conditions, 

sa 2 X 10~^ s, which is much shorter than the accretion timescale tacc of a fluid element, so that weak interactions determine the 
composition. Photons are trapped and are advected with the flow. For neutrinos the situation is more complicated, since the optical 
depth is of order unity. In the outer regions of the disk, they escape freely, whereas for densities larger than « 10^' g cm"^ they 
undergo diffusion on a timescale taa^v w 30 ms. 

2.2. The equation of state and the composition of the fluid. 

For simplicity of presentation, we wiU first assume that all nucleons are free (no a particles are present), and include the necessary 
corrections subsequently. 

Under the conditions described above, the temperature is high enough for electron-positron pair creation. The number of pairs 
thus produced is very sensitive to the degeneracy of the electrons. In fact the number density of pairs is suppressed exponentially 
with increasing degeneracy, because of Fermi blocking. In previous work (Popham, Woosley & Fryer 1999; Narayan, Piran & 
Kumar 2001; Kohri & Mineshige 2002; DiMatteo, Perna & Narayan 2002; Lee, Ramirez-Ruiz & Page 2004) the assumption of full 
degeneracy has been made for simphcity (note also that in some cases the presence of electron-positron pairs has been assumed 
while at the same time retaining a degeneracy pressure term, which is inconsistent). Here we take a different approach, using an 
exact expression for the pressure as a function of the temperature and the chemical potential, valid for arbitrary degeneracy in the 
hmit of relativistic electrons (i.e., p > 10^ g cm"^), due to Blinnikov, Dunina-Barkovskaya & Nadyozhin (1996), namely: 



1 



n-K^hcf 

The number densities of electrons and positrons are related by: 

pYe 1 



(1) 



3.^(..)3-[^^'^^-^^^^>'] 

where is the electron fraction. The chemical potential of species i is denoted by 77, throughout. This expression reduces to the 
well-known Umits when the temperature is low (kT <C r]e, which gives P oc p'^l^, appropriate for a cold relativistic Fermi gas) and 
when it is high (^T » r/^, which gives P oc T'^, when the pressure comes from relativistic electron-positron pairs). The full equation 
of state then reads: 

P = Pad + Pgas + Pe + Pv, (3) 

where 

Prad = ^, (4) 
^gas — ) (5) 

and Pv is the pressure due to neutrinos (discussed below). Here a is the radiation constant, k is Boltzmann's constant and m„ = 
1 .667 X 10"^* g is the atomic mass unit. Since the presence of pairs is automatically taken into account in the expression for Pe, there 
is no alteration to the numerical factor 1 /3 in the expression for P^.^^. For the conditions encountered in the accretion disks presented 
below, gas pressure dominates at the 80% level over the other terms. Regarding the photons, since the temperature is T « 5 MeV, the 
peak in the blackbody spectrum is at ~ 14 MeV. On the other hand, the plasma frequency, Wp, corresponds to Tp « 0.5 MeV. Thus 
the standard expression for radiation pressure may be considered accurate, and plasma effects negligible as far as the photons are 
concerned. 
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The computation of the electron fraction and the chemical potential of electrons follows from the assumption of /^-equilibrium 
between neutrons, protons and electrons, and the condition of charge neutrality. As noted by (Beloborodov 2003), a distinction 
needs to be made to determine the equiUbrium composition depending on the optical depth of the material (the determination of the 
opacities is presented in § 2.4). If it is transparent to its own neutrino emission, we fix equilibrium by equating the capture rates of 
electrons and positrons onto protons and neutrons respectively. For mild degeneracy, as is the case here, this leads to the following 
expression for the electron fraction as a function of the temperature and the electron chemical potential (Beloborodov 2003): 

F. = i+0.487f^Z^y (6) 



2 \ kT _ 

where Q = inin - mpY" ^ 1 .29 MeV. If, however, the material is opaque, and the neutrinos are allowed to diffuse out on a timescale 
shorter than the accretion timescale, then we may write 

Ve + Vp = Vn, (7) 

for the equilibrium composition'*. Since the nucleons are not degenerate, we may use Maxwell-Boltzmann statistics to describe their 
distribution function, and obtain 

^=exp[(e-7?,)/A:r] (8) 
n„ 

for the ratio of proton to neutron number densities. Further, with Ye = np/ (Hp + n„) we arrive at 

1-K 



J, exp[(?7,-0AT]. (9) 

We now have the three equations (2), (3) and (6) or (9) for the three functions T,Ye,T]e and so the system is closed. To allow 
for a transition from the optically thin to optically thick regime, in fact we solve these in a combined form, weighted by factors 
/(Ty) = exp[-ri^] or giTj,) = (1 - exp[-Ti,]). As a practical matter, the internal energy per unit mass, u, is used instead of the pressure 
to solve equation (3), since its variation in time is what is determined in the code, using the First Law of Thermodynamics. We 
finally include the effects of incomplete photodisintegration of a particles by using nuclear statistical equilibrium for the three 
species in,p, a) to fix the mass fraction of free nucleons as (Qian & Woosley 1996): 

^nuc = 22.4(^j (to^J exp(-8.2[10-K]/r). (10) 

Whenever this expression results in X^uc > 1 we set Xnuc = 1 . The corresponding alterations in the previous derivation are simple, and 
lead to the full set of equations which we write for each gas element and at each time. Note that these do not allow for an explicit 
solution, so an iterative scheme is used in each instance. The internal energy per unit mass is 

w = 3 + — , (11) 

p 4 2w„ 



/3-equilibrium gives 

Ye = (1 ~-Xi[uc)/2+-Xni 

and charge neutraUty impUes 



\0.487 ^2/^-^^ 



2 I kT 

pYe _ 1 



/(T.) + [l+exp({7?,-e})/A:r)k(T.) , (12) 



[r]l + r]eTT\kTf], (13) 



as given by equation (2). The modifications to the equation for f3 equilibrium simply reflect the fact that if the fluid is composed 
primarily of a particles, there will be one neutron per proton, and hence Y^^ 1/2. 

2.3. Neutrino emission and photodisintegration losses 

Several processes contribute to the emission of neutrinos. First, we take into account electron and positron capture onto free nu- 
cleons using the tables of Langanke & Martinez-Pinedo (2001). This is an improvement over calculations done previously (Popham, 
Woosley & Fryer 1999; Narayan, Piran & Kumar 2001; Kohri & Mineshige 2002; DiMatteo, Perna & Narayan 2002; Lee, Ramirez- 
Ruiz & Page 2004), where assumptions concerning the degeneracy were made. The rates are thus accurate over the entire disk, 
whether the degeneracy is significant or not. The tables cover the following ranges: 1 < logipYglg cm"-']) < 11; 7 < logr[K] < 11. 
A bilinear interpolation in the log pF^, - log T plane is performed using the table to obtain the cooling rate for a given mass element. 
The result is then multiplied by the mass fraction of free nucleons, Xn^c, since we are not considering capture of electrons and protons 
by a particles. Second, the annihilation of electron-positron pairs is a source of thermal neutrinos, and the corresponding cooling 
rate is computed using the fitting functions of Itoh et al. (1996). These cover the range 9 < logiplg cm~^]) < 12 in density and 
10 < log ^[K] < 1 1 in temperature. Finally, nucleon-nucleon bremsstrahlung, qff, and plasmon decay, ^'piasmon. are considered, with 
rates given by: 

(iff = 1 .5 X 10^^ r/i^ pj^ erg s"' cm'^ , (14) 

This expression neglects the neutrino chemical potential, r),^, and is strictly valid oiHy when there are equal numbers of neutrinos and anti-neutrinos. Otherwise the 
electron pressure, bulk viscosity and Joule-Uke heating due to the induced lepton current will be affected (Burrows, Mazurek & Lattimer 1981; Socrates et al. 2004). 
We neglect all these corrections in the present treatment. 
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(Hannestad & Raffelt 1998) and 

^piasmon = 1-5 X 10'%\ 7^ exp(-7p)(l +7^) ^+ erg s-i cm-3, (15) 

where 7^ = 5.5 x lO~^^yi^^'^ + 3[r]e/kT]'^)/3 (Ruffert, Janka & Schafer 1996). For the conditions encountered in the disk, capture by 
nucleons completely dominates all other processes, and is the main source of cooling. Finally, the creation and disintegration of a 
particles leads to a cooling term in the energy equation given by: 

^phot = 6.8 X lO^** erg s"' cm'^. (16) 
2.4. Neutrino opacities 

The material in the disk is dense enough that photons are completely trapped, and advected with the flow. For neutrinos however, 
the situation is more complicated, since the outer regions are transparent to them, while the inner portions are opaque. 

The main source of opacity is scattering off free nucleons and a particles, with cross-sections given by (Tubbs & Schramm 1975; 
Shapiro & Teukolsky 1983) 

and ^ 

cT^=ao(^] [(4sin2^„)]', (18) 

where cto = 1.76 x 10"^ cm^ and 9^ is the Weinberg angle. Since capture processes dominate the neutrino luminosity, we as- 
sume that the mean energy of the neutrinos is roughly equal to the Fermi energy of the partially degenerate electrons » |T/g = 

43(FePi2)^/^ MeV, where pn = p/Wh 

cm ^. So the mean free path is /;/ = l/(nNCN+Wafa), andna being the number densities 
of free nucleons and a particles respectively. To compute an optical depth, we take Ti, = H/l,y, with H a typical scale height of 
the disk. Inspection of the disk shape (through the density contours in the inner regions) and an estimation of the scale height with 
~ p/Vp suggests that H ~ Kr, with k a constant of order unity, so that our final, simplified expression for the optical depth reads 
Ti, = Kr/li,. With the above expressions for the cross section as a function of density and electron fraction, and defining rj = r/lO^cm 
we arrive at 

= 186.5 KPif y2/3 [Xnuc + 3.31(1 -X„ue)/4] . (19) 

The term which depends on Xnuc reflects the two-species composition we have assumed. Its influence is Umited, since when Xnuc 
varies from to 1, the optical depth is altered by a factor 1.2, all other values being equal. In practice, and for the compositions found 
in the disks, this expression corresponds to having an optical depth of unity at approximately 10"g cm"^, which can be thought of 
as a neutrino-surface, where the last scattering occurs and neutrinos leave the system. Accordingly, we modify the expressions for 
the neutrino emission rates, suppressing it in the opaque regions, and enhancing the pressure through 

/ du\ ( du\ 

and 

P, = ^ar4[l-exp(-r,)], (21) 

where {du / dt){) = X] qi/ p is the unmodified energy loss rate (in erg g~' s~') calculated from the rates given in § 2.3. 
The total neutrino luminosity (in erg s"') is then computed according to 

I ' [^•^•^"'■^P'^™°""'"^P'^"'"^'=^p] exp(-r,.)Jm. (22) 

2.5. Numerical method and initial conditions. 

For the actual hydrodynamical evolution calculations, we use the same code as in previous work (Lee & Ramirez-Ruiz 2002), and 
refer the reader to that paper for the details. This is a two dimensional (cyhndrical coordinates {r,z) in azimuthal symmetry) Smooth 
Particle Hydrodynamics code (Monaghan 1992), modified from our own 3D version used for compact binary mergers (Lee 2001). 
The accretion disk sits in the potential weU of a black hole of mass Mbh, and which produces a Newtonian potential $ = -GMBn/r. 
Accretion is modeled by an absorbing boundary at the Schwarzschild radius rg^h = '^GM-^^/c^, and the mass of the hole is updated 
continuously. The transport of angular momentum is modeled with an a viscosity prescription, including all components of the 
viscous stress tensor (not only trtf,). The self-gravity of the disk is neglected. 

The main modifications from the previous version are in the implementation of a new equation of state §2.2, neutrino emission, 
§2.3 and the treatment of neutrino optical depths, §2.4. 

As done previously (Lee & Ramirez-Ruiz 2002; Lee, Ramirez-Ruiz & Page 2004), the initial conditions are taken from the final 
configuration of 3D calculations of black hole-neutron star mergers, after the accretion disk that forms through tidal disruption of 
the neutron star has become fairly symmetric with respect to the azimuthal coordinate. We show in Table 1 the parameters used in 
each of the dynamical runs included in this paper. 
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3. RESULTS 

As the dynamical evolution calculations begins, the disks experience an initial transient, which is essentially numerical in origin, 
and is due to the fact that the configuration is not in strict hydrostatic equilibrium (recall that it is obtained from azimuthally averaging 
the results of 3D calculations). This transient, and its effects, are negUgible, and it is essentially over in a hydrodynamical timescale 
(w 2 ms). After that, the disk proceeds to evolve on a much longer timescale, determined by accretion onto the central black hole. We 
first present the details of the spatial structure of the disk due to fundamental physical effects, and then proceed to show the temporal 
evolution and associated transitions. 

3.1. Disk structure 

The fundamental variable affecting the instantaneous spatial structure of the disk is the optical depth to neutrinos, r,y, since it 
determines whether the fluid cools efficiently or not. All quantities show a radical change in behavior as the threshold Tj^ = 1 is 
passed. Figure 1 shows color-coded contours in a meridional sUce of the thermodynamical variables, while Figure 2 displays the 
run of density and entropy per baryon along the equator, z = 0, where the density and temperature are highest (we will refer here 
mainly to run a2M, unless noted otherwise). The density increases as one approaches the black hole, varying as p oc r"^ in the outer 
regions, where <C 1. The critical density for opaqueness is reached at ~ 10^ cm, and for r < r*, p cx Since the energy that 
would otherwise be lost via neutrino emission remains in the disk if the cooling is suppressed by scattering, the density does not rise 
as fast. In fact, one can note this change also by inspecting the scale height H ~ P/VP, which scales as ///r oc r for r > and 
H/r oc const for r < (and thus in terms of surface density the change is from S oc to S oc const at r*). The opaque region of 
the disk remains inflated to a certain extent, trapping the internal energy it holds and releasing it only on a diffusion timescale. At 
densities p ~ 10^^ g cm"^ the optical depth can reach t,^ ~ 100, and thus the suppression of cooling is dramatic, as can be seen in 
Figure 1, where the contours of and q are shown. The flattening of the entropy profile is related to the change in composition in 
the optically thick region and the occurrence of convection (see § 3.2 below). 

As the density rises in the inner regions of the disk, the electron fraction initially drops as neutronization becomes more 
important, with the equilibrium composition being determined by the equality of electron and positron captures onto free neutrons 
and protons (see equation (6) and the discussion preceding it). The lowest value is reached at r ~ r*, where ~ 0.03. Thereafter 
it rises again, reaching ~ 0.1 close to the horizon. Thus flows that are optically thin everywhere will reach a higher degree of 
neutronization close to the black hole than those which experience a transition to the opaque regime. The numerical values for the 
electron fraction at the transition radius and at the inner boundary are largely insensitive to a, as long as the transition does occur. 

The baryons in the disk are essentially in the form of free neutrons and protons, except at very large radii and low densities 
(r > 5 X 10^ cm and p <5x lO^g cm~^) where a particles form. Figure 3 shows the region in the density-temperature plane where 
the fluid lies, for runs a2M and alM at two different times (and color coded according to the electron fraction). Most of the gas 
lies close to the line determining the formation of Helium nuclei (see equation 10), but does not cross over to lower temperatures. 
The reason for this is that the energy that would be released by the creation of one Helium nucleus (28.3 MeV) would not leave the 
disk (recall that we consider neutrino emission as the only source of cooling) and thus immediately lead to the photodisintegration 
of another a particle. An equiUbrium is thus maintained in which the gas is close to HeUum synthesis, but this does not occur. In 
the opaque regime, the gas moves substantially farther from the transition Une to HeUum, since it cannot cool efficiently and, as 
mentioned before, the density does not rise as quickly. Figure 3 also shows clearly why one must make use of an equation of state 
which takes into account properly the effects of arbitrary degeneracy of the electrons and positrons. The solid straight line marks 

1/3 

the degeneracy temperature as a function of density, given by kT = 1-1 Pi[ MeV. The disk straddles this line, with a degeneracy 
parameter rje/kT ~ 2-4 in the inner regions, and rje/kT Ri 1 at lower densities. Thus making the approximation that the electrons in 
the flow are fully degenerate is not accurate. 

Modeling the accretion flow in two dimensions (r,z), without the assumption of equatorial symmetry allows one to solve clearly 
for the vertical motions in the disk, something which is not possible when considering vertically-integrated flows. It was pointed 
out by Urpin (1984) that a standard a viscosity could lead to meridional flows in which Vr changed sign as a function of height 
above the mid plane, leading to inflows as well as outflows. Kita (1995), Kluzniak & Kita (2000) and Regev & Gitelman (2002) 
later considered a similar situation, and found that for a range of values in a, the gas flowed inward along the surface of the disk, 
and outward in the equatorial regions. In previous work (Lee & Ramirez-Ruiz 2002), we found this solution in disk flow, with large 
scale circulations exhibiting inflows and outflows for high viscosities (a ~ 0.1), and small-scale eddies at lower values (a < 0.01). 
How vertical motions affect the stability of accretion disks and may lead to the transport of angular momentum is a question that has 
become of relevance in this context (Arlt & Urpin 2004). Here we show in Figure 4 the velocity field and magnitude of meridional 
velocity for run a2M. The small-scale eddies are clearly visible, with their strength usually diminishing as the disk is drained of 
matter and the density drops. The effect of a different value of a can be seen in the top row of Figure 5, where the comparison 
between a = 0.1 anda = 0.01 is made. The large, coherent Unes of flow aimed directly at the origin for a = 0.1 correspond to inflow, 
while the lighter shades along the equator at larger radii show outflowing gas. 

The instantaneous structure of the disk concerning the density, temperature and composition profiles is largely independent of the 
viscosity, as long as there is enough mass to produce the optically thin/optically thick transition. Even in the optically thin regime, 
the structure cannot be determined analytically following the standard arguments applied to thin, cool disks of the Shakura-Sunyaev 
type. The reason for this is the following: a central assumption in the standard solution is that the disk is cool, i.e., kT /nip-^ GM^n/r. 
This comes from the requirement that all the energy dissipated by viscosity, qo, be radiated away efficiently, and produce the observed 
flux, F . In our case, the disk is most certainly not cool, as it originates from the tidal disruption of a neutron star by a black hole (or 
possibly the coalescence of two neutron stars, and subsequent collapse of the central mass to a black hole). The gas that constitutes 
the disk is dynamically hot because it was in hydrostatic equilibrium in a self-gravitating configuration, where U ~ —W and has 
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not been able to release this internal energy (the merger process itself may lead to additional heating). A further deviation from the 
standard solution is that the pressure support in the disk leads to a rotation curve that is slightly sub-Keplerian. The dissipation by 
viscosity is in fact smaller than the cooling rate over much of the disk, and the released luminosity comes from a combination of 
viscous dissipation and the store of internal energy given to the disk at its conception. 

3.2. Lepton-driven convection 

The classical requirement for convective instabiUty in the presence of entropy and composition gradients, as well as rotation, is 
the Solberg-Hoiland criterion (Tassoul 1978): 

N^ + Loj<0, (23) 
where = 4Q^ + rdff'/dr is the radial epicyclic frequency (fl being the angular velocity), and 
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1 /dP\ ds \ f dP\ dYe 
P ydsjy dr'^P\dYe)^~dr 



(24) 



Ye ' 

is the Brunt-Vaisala frequency (see Lattimer & Mazurek 1981; Thompson, Quataert & Burrows 2005). The adiabatic index 7 = 
dlnP /dlnp w 5/3 in our case, since Pgas is the most important contribution to the total pressure. Thus a region may be convectively 
unstable because of a composition gradient, or an entropy gradient, or both. Strictly speaking, here one should consider the total 
lepton fraction F;. We do not consider the transport of neutrinos in detail, as already mentioned (see § 2.4), and do not calculate 
explicitly the contribution of neutrinos to this 7;. For what follows, we will thus simply take Yg to represent the behavior of the full 
Yi. The origin of convection in the lepton inversion zone can be understood as follows (Epstein 1979). Consider a fluid element in 
the lepton inversion zone that is displaced in the outward direction and then comes to pressure equilibrium with its surroundings. 
The displaced element, which is lepton rich relative to its new surroundings, attains the ambient pressure at a lower density than the 
surrounding fluid, because the pressure depends directly on the lepton number, and thus tends to drift outwards. By the same token, 
an inwardly displaced fluid element in the lepton inversion zone is depleted in leptons relative to its new surroundings and thus tends 
to sink ^. When temperature gradients are allowed for, the displaced fluids, in general, have temperatures which differ from those of 
their surroundings. These variations tend to promote stability or instabihty depending on whether the existing temperature gradient 
is less than or greater than the adiabatic gradient, respectively. 

In the present scenario then, the transition to the optically thick regime leads to instability, because initially ds/dr < 0, and also 
dYe/dr < (see also Figure 2). The entropy profile is then flattened by efficient convective mixing, and the lepton gradient inversion 
is due to the different way in which the composition is determined through weak interactions once the neutrinos become trapped. 
In a dynamical situation such as the one treated here, convection tends to erase the gradients which give rise to it. Accordingly, 
the sum of the corresponding terms in equation (23) tends to zero once the simulation has progressed and convection has become 
established in the inner disk (note that in our case the rotation curve there is sub-Keplerian, with O oc r~^l^, so that Q. and lj,- are not 
equal). This is analogous to what occurs in proto neutron stars following core collapse (Epstein 1979; Burrows & Lattimer 1986; 
Thompson & Duncan 1993) and has actually been confirmed in numerical simulations of such systems (Janka & Muller 1996), where 
the convection leads to strong mixing. For a neutrino-driven convective luminosity Leon at radius R and density p, the convective 
velocity may be estimated by standard mixing length theory as (see e.g., Thompson & Duncan 1993) 

v,o„ ^ 3.3 XIOM T7uA-T Utt^ cms-', (25) 



10^2 erg s-iy \^10'2gcm-V V20km 
where we have used the fact that gas pressure dominates in the fluid. The assumption of spherical symmetry imphcit in this expression 
is not strictly met in our case, but it may provide us nevertheless with a useful guide. The overturn time of a convective cell is then 
given by fcon ^ Ip/vcon ^10 ms, where Ip ^ 20 km is the mixing length, usually set to a pressure scale height. In our calculations, 
we see that the magnitude of the meridional velocity | v| = ^yvj + v^ decreases as r decreases, then rises again as the opaque region is 
reached, reaching |v| ~ 10^ cm s~\ at r ~ 20 km. The associated turnover times are thus tcon — lp/\v\ — 20 ms, in good agreement 
with the estimate made above. 

To isolate this effect from that due to viscosity (which generates the meridional circulations mentioned previously) we have 
performed a simulation in which a = (run alM in Table 1). This calculation shows essentially no accretion onto the black hole 
except for a small amount of gas transferred at early times (because the initial condition is not in strict equilibrium). The result 
concerning the profile of Ye as a function of radius is as we have described above, i.e., increasing neutronization as r decreases, until 
the opaque region is reached, followed by an increase in Ye as the black hole boundary is approached. The revealing difference lies 
in the time evolution of this profile. As mentioned above, convection generates motions which tend to eliminate the composition 
gradient which drives it. In the presence of viscosity, matter is continuously transported radially, and the gradient is not erased 
entirely. Height integrated profiles of Ye{r) at 50 ms and 200 ms show a similar behavior at small radii. When viscosity is removed, 
the initial composition gradient is gradually softened until it disappears in the innermost regions (see Figure 7). The disk then has 
a nearly constant electron fraction along the equator for r < , with a sharp transition region leading to an increase with radius in 
the transparent regime. This may resemble the convection dominated accretion solution found analytically by Quataert & Gruzinov 
(2000) and numerically by several groups (Stone, Pringle & Begelman 1999; Narayan, Igumenshchev & Abramowicz 1999), where 
convection transports angular momentum inward, energy outward, and gives a radial profile in density a We find a different 
power law, with p oc r"'. There are several factors which may account for this difference: the disk vs. spherical geometry; the initial 
condition with a large amount of internal energy; and the limited, but present cooUng rate on the boundaries of the flow. 

' This is similar to the thermosolutal convection which would occur in a nonnal star if there were an inwardly decreasing molecular weight gradient or composition 
inversion (Spiegel 1972). 
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3.3. Neutrino diffusion effects 

Since the neutrinos are diffusing out of the fluid (the mean free path is small compared to the size of the system in the optically 
thick regime), one would expect a corresponding viscosity through Q"^ (see, e.g. Burrows, Mazurek & Lattimer 1981). This needs 
to be compared to the corresponding viscosity generated by the assumed a prescription for consistency. We may assume C"** = 
{\/3)Uvlv/c, where Uv ~ aT^ is the energy density associated with neutrinos, and is their mean free path (see § 2.4). Thus 

C™«2x 10^°gcm-i s-i, (26) 

for conditions near the equatorial plane, z = 0. This value will increase in the outer regions, since the mean free path becomes larger 
as the density drops. At the neutrino-surface, where « 1, the mean free path is 1^ ~ 10^ cm, and so an averaged value of the 
viscosity over the entire neutrino-opaque region will result in an increase of about one order of magnitude over the estimate given 
in equation (26) (note however, that strictly speaking the diffusion approximation is no longer valid in the outer regions, so this must 
be interpreted with care). For the a prescription, C.^ = pctc^l^- The rotation curve is not too far from Keplerian, and scaling this 
expression to typical values we find 

Ca: = 1 .3 X 10^2 c^g a-2 r^M'^'^ p^g cm"' s"' . (27) 

The effects on angular momentum transport are thus a fuU two orders of magnitude below those arising from our viscosity prescription 
(for a = 10"^). To put it another way, a lower limit for the viscosity under these conditions (and in the optically thick portion of the 
disk) would be a ^ 10"^. An alternative analysis would be to consider the corresponding timescales induced by this viscosity. Since 
fvis oc R{p/ Q^, smaller viscosities imply longer timescales, as expected. 

3.4. Stability 

Aside from convection, several general criteria for stability can be analyzed in our case. Evidently, since we are performing 
dynamical calculations, any instability that arises wiU quickly lead to a change in structure. It is nevertheless instructive to consider 
the corresponding conditions within the disk. We first consider the Toomre criterion, comparing gravitational and internal energies, 
with 

Qt = ^-^, (28) 

where ujr is the local epicyclic frequency (essentially equal in this case to the angular frequency O), c, is the local sound speed and 
S is the surface density. We find gx > 1 in all cases throughout the calculations, and thus that the disks are stable in this respect 
(Figure 8a shows a typical profile of 2x1 ''I)- 

Previous studies (Narayan, Piran & Kumar 2001; Kohri & Mineshige 2002) have shown that neutrino cooled disks are thermally 
unstable if radiation pressure dominates in the flow, and stable otherwise, although the effects of neutrino opacities were not consid- 
ered. The criterion for stability can be written in this case as 

' ^ ^ < , (29) 



\ dhiT J ^ \ dlnT 

where and Q~ are the volume heating and cooling rates. Essentially, in order to be stable the disk must be able to get rid of any 
excess internal energy generated by an increase in the heating rate. We show in Figure 8b typical values for and Q~ as functions 
of the central (equatorial) temperature, for run a2M at f = 50 ms. In the optically thin region the disks are thermally stable, because 
the pressure is dominated by the contribution from free nucleons and cooling is efficient. In the optically thick region the cooling is 
greatly suppressed and the criterion would indicate that the disk becomes thermally unstable. This is reasonable, since with optical 
depths Tj^ « 10- 100, essentially no direct cooling takes place, and dissipation is not suppressed. This is why the disk is geometrically 
thick in the inner regions, as already described above (§3.1). The balance that gives a quasi-steady state structure is achieved mainly 
through the diffusion of neutrinos, since the diffusion timescale is shorter than the accretion timescale. 



3.5. Disk evolution 

The evolution of the disk on long timescales is determined by the balance between two competing effects: on one hand, viscosity 
transports angular momentum outwards, matter accretes and the disk drains into the black hole on an accretion timescale ?acc- The 
trend in this respect is towards lower densities and temperatures. On the other hand, cooHng reduces pressure support and leads to 
vertical compression, increasing the density. The internal energy of the fluid is released on a cooling timescale, fcooi- The presence 
of an optically thick region in the center of the disk limits the luminosity (dominated by electron and positron capture onto free 
nucleons) to « few x lO^^erg s~S and the initial internal energy is « 10^^ erg, so fcooi — 0.1 s. If face > fcooi, the disk will cool 
before its mass or internal energy reservoir is significantly affected by accretion onto the black hole. The maximum density (shown 
in Figure 9 for runs alM, a2M and a3M) actually increases slightly due to vertical compression, and subsequently drops once mass 
loss through accretion dominates (the ^ 10 ms delay at the start, during which it is approximately constant, is simply the sound 
crossing time across the optically thick region of the disk). The accretion rate onto the black hole, Mbh, the accretion timescale 
tacc = A/flH/A^disk and the total neutrino luminosity L^, are shown in Figures 10, 11 and 12. They aU show the same qualitative 
behavior, remaining fairly constant (or changing slowly) for an accretion timescale, and abruptly switching thereafter. The accretion 
timescales are approximately 0.5 s and 5 s for a = 0.01,0.001 respectively. 

For high viscosity, a = 0.1, the transport of angular momentum is vigorous, and the black hole quickly accretes a substantial 
amount of mass (0.16 Mq within the first 100 ms). The accretion timescale is ?acc — 50 ms. The circulation pattern consists of 
large-scale eddies, with // r. In fact there is essentially one large eddy on each side of the equatorial plane, with mass inflow along 
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the surface of the disk, and an equatorial outflow. Part of the outflowing gas moves away from the equator and reverses direction 
close to the surface of the disk, contributing to the inflow. For an intermediate viscosity, a = 0.01, the intensity of the circulations 
is smaller, but also, the eddies are smaller, with several of them clearly occurring in the disk at once. The transport of angular 
momentum being less vigorous, the accretion rate is substantially smaller than in the previous case. Finally, for a yet lower viscosity, 
a = 0.001, the trend continues, and the eddies become smaller still. In these last two cases, angular momentum transport is so low 
that very high densities are maintained in the central regions of the disk for a large number of dynamical times, contrary to what is 
seen in the high viscosity case. A simple way to quantify how much of the mass flow is actually making it to the central black hole 
is to measure the fraction of mass at any given radius that is moving inwards, |Min|/(|Min| + |Mout|)> where we have divided the mass 
flow rate into two height-integrated parts, one with Vr < and another with Vr > 0. For example, for run a2M at ? = 100 ms and 
r = r* it is approximately 1/3. The ratio tends to unity only in the innermost regions of the disk, and shows the true black hole mass 
accretion rate, plotted in Figure 10 for the same cases. 

We now turn our attention to the neutrinos, which are the main source of cooUng. The results in this case are markedly different 
from what we initially found (Lee & Ramirez-Ruiz 2002), simply because of the new equation of state, the more realistic cooling 
rates and the approximate computation of opacities. They are in general agreement with the preliminary results we presented before 
(Lee, Ramirez-Ruiz & Page 2004), which used a less detailed equation of state than the one shown here. 

The optically thick region is present in every disk at the start of the calculation. As already mentioned, this has two important 
effects: enhancing the pressure and suppressing the neutrino emission. This reflects upon the total luminosity, since the suppression 
occurs precisely in the hottest regions, where most of the energy would otherwise be released. The most important qualitative 
difference between the runs presented here is that for a high viscosity (a = 0.1, runs alM and aim), the disk is drained of mass so fast 
that it has no chance to cool (4cc < fcooi) and release most of its internal energy. It is in fact advected into the black hole. Moreover, 
the optically thick region disappears entirely from the disk (see Figure 5) by ? = 40 ms. The emission is then no longer suppressed 
(see the contours of q in the same figure), and the disk radiates at the maximum possible cooling rate. The thermal energy content 
of the disk is so large, that as it thins, the luminosity actually increases briefly around t = 10 ms before dropping again. The drop at 
late times follows an approximate power law, with oc . In this interval the energy release comes from a combination of residual 
internal energy and viscous dissipation within the disk. For lower values of the viscosity, a < 0.01, the central regions of the disk 
remain optically thick throughout the calculations. Note that reducing the disk mass by a factor of five (as was done for runs aim, 
a2m, a3m) does not affect these overall conclusions. The fluid is simply compressed into a smaller volume, and thus the densities 
and temperatures that are reached are similar than for the high mass runs. For run a3M the neutrino luminosity is practically constant 
at 3 X 10^^ erg s"' for / > 100 ms. Table 2 summarizes the typical disk mass, energy density, accretion rate, luminosity, and duration 
and energetics of neutrino emission for aU runs. 

4. SUMMARY, CONCLUSIONS AND ASTROPHYSICAL IMPLICATIONS 

4.1. Summary 

We have performed two-dimensional hydrodynamical simulations of accretion disks in the regime of hypercritical accretion, where 
neutrino emission is the main cooUng agent. The disks are assumed to be present around a stellar-mass black hole, and are evolved for 
a few hundred dynamical timescales. We have paid particular attention to the relevant microphysical processes under the conditions 
at hand, and used a detailed equation of state which includes an ideal gas of a particles and free baryons in nuclear statistical 
equilibrium and a relativistic Fermi gas of arbitrary degeneracy. The composition of the fluid is determined by weak interactions. 
The density and temperature are such that the iimer regions of the disk become opaque to neutrinos, and this is taken into account in 
a simple approximation. 

4.2. Conclusions 

Our main conclusions can be summarized as follows: 

• Once the fluid becomes photodisintegrated into free nucleons, neutronization becomes important and lowers the electron 
fraction substantially below 1/2, with the electron fraction Ye reaching w 0.05 at its minimum. This value, however, does not 
occur in the immediate vicinity of the black hole, but rather at the transition radius where the fluid becomes optically thick 
to its own neutrino emission, w 10^ cm. At smaUer radii, the electron fraction rises again, reaching w 0.1 close to the 
horizon. 

• Neutrino trapping produces a change in composition and an inversion in the electron fraction. The associated negative 
gradient in Ye{r) induces convective motions in the optically thick region of the disk. This is analogous to what presumably 
occurs following core collapse, in a proto-neutron star and its surrounding envelope. Due to the radial flows induced by 
viscosity, convection is unable to suppress this composition gradient. It would appear, however, that the entropy per baryon 
in the opticaUy thick region is very close to being constant, with s/k»6. 

• The spatial structure of the disk is characterized by the transition radius r, where the material becomes optically thick. 
For r > the disk cools efficiently, whereas for r < the emission is suppressed and the fluid is unable to cool directly 
(although it does so on a neutrino diffusion timescale). This leads to larger pressures and a more moderate rise in density. 

• The temporal evolution of the disk is determined by the balance between accretion and neutrino emission. For low viscosities 
(a < 0.01), the disk is able to cool in a quasi steady state and radiate its internal energy reservoir. This lasts for approximately 
0.1-0.4 s, with w 10^^ erg s~^ Thereafter the typical luminosity and density quickly decay. For large viscosities (a ~ 0.1) 
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the disk is drained of mass on an accretion timescale which is shorter than the cooUng timescale, and the internal energy of 
the disk is essentially advected into the black hole. An interesting result in this case is that as the disk becomes transparent 
before being engulfed by the hole, it undergoes a re-brightening, as some of the stored energy escapes. 

• The total energy output in neutrinos is Ei, ss 10^^ erg, over a timescale of « 200 ms. The typical accretion rates are 
PS O.IM0 s~\ and neutrino energies are « 8 MeV at r*. Energy densities in the inner regions of the disk are w lO^^erg cm~^. 

4.3. Discussion 

There are two main ingredients in the results presented here that contrast with those available previously in the literature concerning 
the steady state structure of neutrino cooled accretion disks. The first is the ability to dynamically model the evolution of the 
system for hundreds of dynamical timescales, taking the previous history of the fluid into consideration through the choice of initial 
conditions. This allows us to consider the energetics on more relevant timescales (cooling, viscous) than the dynamical one accessible 
in three dimensional studies. The second is the realization that the structure of the disks is affected qualitatively by the presence of 
an optically thick region at high densities. In this sense the situation is similar to that encountered following massive core collapse, 
where convection occurs (recent work assessing the relative importance of the MRI, neutrino and convection driven viscosity in 
rotating collapsing cores has been reported by Akiyama et al. 2003; Thompson, Quataert & Burrows 2005, we comment further on 
this issue in §4.3.2). 

Clearly there is room for improvement in the results presented here. To begin with, our expressions for the effect of a finite optical 
depth on cooling and pressure are too simple, and attempt only to capture the essential physical behavior of the system. There is an 
important region of the disk where the optical depth is not large enough to consider the diffusion approximation, and where more 
detailed transport effects ought to be considered. We have not separated the neutrino variables into three species, which would be 
more rigorous (e.g., Yokosawa et al. 2004, have noted that in such flows, the electron neutrinos Vg might be preferentially absorbed 
with respect to electron anti-neutrinos Ve, thus affecting the neutrino annihilation luminosity above the surface of the disk in an 
important way). We have only considered the effects of coherent scattering off nucleons and a particles for opacity purposes. This 
specifically ignores absorptive scattering, which would: (i) lead to a modification of the neutrino emergent spectrum; (ii) produce 
heating of the fluid. The latter may be particularly important concerning the driving of powerful winds off the surface of the disk, 
and needs to be addressed more carefully in the context of GRBs. As in previous work, we have chosen to maintain the Newtonian 
expression for the potential well of the central black hole for ease of comparison to previous results (this will be addressed in future 
work). 

4.3.1. Gamma-Ray Bursts 

The sudden release of gravitational binding energy of a neutron star is easily sufficient to power a GRB. The minimum energy 
requirement is 10^'(ArJ/47r) erg if the burst is beamed into a solid angle Afi. As discussed in §1, a familiar possibiUty, especially for 
bursts belonging to the short duration category, is the merger of a neutron star - black hole, or double neutron star binary by emission 
of gravity waves, which, as illustrated here, is likely to generate a black hole surrounded by a lower mass accretion disk. How is the 
available rotational and gravitational energy converted into an outflowing relativistic plasma? A straightforward way is that some of 
the energy released as thermal neutrinos is reconverted, via collisions outside the disk, into electron-positron pairs or photons. The 
neutrino luminosity emitted when disk material accretes via viscous (or magnetic) torques on a timescale Af ~ 1 s is roughly 

^.~2xlO«(ilSfV^)"'e,gs-. (30, 



.O.IMqJ V Is/ 

for a canonical radiation efficiency of 0.1. During this time, the rate of mass supply to the central black hole is of course much 
greater than the Eddington rate. Although the gas photon opacities are large, the disk becomes sufficiently dense and hot to cool 
via neutrino emission. There is in principle no difficulty in dissipating the disk internal energy, but the problem is in allowing these 
neutrinos to escape from the inflowing gas. At sufficiently low accretion rates, a < 0.01, we find that the energy released by viscous 
dissipation is almost completely radiated away on a timescale given by tcooX ~ £int/^i^ ~ 0.1 s. In contrast, for a higher mass supply, 
a > . 1 , energy advection remains important until the entire disk becomes optically thin. The restriction on the cooling rate imposed 
by high optical depths is key because it aUows the energy loss to be spread over an extended period of time during which the neutrino 
luminosity stays roughly constant. This gives a characteristic timescale for energy extraction and may be essential for determining 
the duration of neutrino-driven short GRBs (Lee, Ramirez-Ruiz & Page 2004). 

Neutrinos could give rise to a relativistic pair-dominated wind if they converted into pairs in a region of low baryon density (e.g. 
along the rotation axis, away from the equatorial plane of the disk). The vD e"^e~ process can tap the thermal energy of the torus 
produced by viscous dissipation. For this mechanism to be efficient, the neutrinos must escape before being advected into the hole; 
on the other hand, the efficiency of conversion into pairs (which scales with the square of the neutrino density) is low if the neutrino 
production is too gradual. Typical estimates suggest a lower bound^ of L^^p ^ IO'-'L,^ (e.g Rosswog & Ramirez-Ruiz 2003). If the 
pair-dominated plasma were coUimated into a solid angle Aft then of course the apparent isotropized energy would be larger by a 
factor (47r/AO), but unless Aft is < lO"'^ this may fail to satisfy the apparent isotropized energy of 10^^ ergs implied by a redshift 
z = 1 for short GRBs. 

One attractive mechanism for extracting energy that could circumvent the above efficiency problem is a relativistic magneto 
hydrodynamic (MHD) wind (Usov 1992; Thompson 1994). Such a wind carries both bulk kinetic energy and ordered Poynting 

* This estimate, however, assumes that the entire surface area emits close to a single temperature black-body. It should be noted that if the dissipation takes place in 
a corona-Uke environment, the efficiency may be significantly larger (Ramirez-Ruiz & Socrates 2005). 
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flux, and it possible that gamma-ray production occurs mainly at large distances from the source (Duncan & Thompson 1992; 
Usov 1994; Thompson 1994). A rapidly rotating neutron star (or disk) releases energy via magnetic torques at the rate L^ag ~ 
10*'^ B^^PZ^Rl erg s"', where P = IQ-^P. 3 s is the spin period, and B = 10'^B|5 G is the strength of the poloidal field at a radius 
R = lQ^R(i cm. The last stable orbit for a Schwarzschild hole lies at a coordinate distance R = 6GM/c^ = 9(M/Mq) km, to be 
compared with R = GM/c^ = 3/2(M/M©) km for an extremal Kerr hole. Thus the massive neutron disk surrounding a Schwarzschild 
black hole of approximately 2Mq should emit a spin-down luminosity comparable to that emitted by a milUsecond neutron star. 
A similar MHD outflow would result if angular momentum were extracted from a central Kerr hole via electromagnetic torques 
(Blandford & Znajek 1977). The field required to produce L^ag > 10''' erg s"' is colossal, and may be provided by a helical dynamo 
operating in hot, convective nuclear matter with a millisecond period (Duncan & Thompson 1992). A dipole field of the order of 
lO'^ G appears weak compared to the strongest field that can in principle be generated by differential rotation (^ 10'^[f /I ms]~' G), 
or by convection (~ 10'* G), although how this may come about in detail is not resolved. We examine in more detail the possible 
generation of strong magnetic fields below in §4.3.2. 

Computer simulations of compact object mergers and black hole formation can address the fate of the bulk of the matter, but there 
are some key questions that they cannot yet tackle. In particular, high resolution of the outer layers is required because even a tiny 
mass fraction of baryons loading down the outflow severely limits the attainable Lorentz factor - for instance a Poynting flux of 
10^^ erg could not accelerate an outflow to F > 100 if it had to drag more than ^ 10~^ solar masses of baryons with it. One further 
effect renders the computational task of simulating jet formation even more challenging. This stems from the likelihood that the high 
neutrino fluxes ablate baryonic material from the surface of the disk at a rate (Qian & Woosley 1996) 

5/3 



M,^5xlO-^( ^Q32e;gs-i ) (31^ 



A rest mass flux hmits the bulk Lorentz factor of the wind to 



r = ^""g = \q( ^ ( ^ 1 (32) 

^ H4 ^2 \ 1052 



Mr/C^ V 10^^ erg s-^ J \5 x IO-^MqS" 
If one assumes that the external poloidal field strength is limited by the vigor of the convective motions, then the spin-down luminosity 
scales with neutrino flux as Lniag oc oc vl^^ oc lI!^, where Vcon is the convective velocity. The ablation rate given in equation (31) 
then indicates that the limiting bulk Lorentz factor of the wind decreases as L^' . Thus the burst luminosity emitted by a magnetized 
neutrino cooled disk may be self-limiting. The mass loss would, however, be suppressed if the relativistic wind were coUimated into 
a jet. This suggests that centrifugally driven mass loss will be heaviest in the outer parts of the disk, and that a detectable burst may 
be emitted only within a certain solid angle centered on the rotation axis (see e.g., Rosswog & Ramirez-Ruiz 2003). 

4.3.2. Generation of Strong Magnetic Fields 

There is also the question of magnetic fields, which we have not included, but should obviously be considered. The field in a 
standard disk is probably responsible for viscous stresses and dissipation, through the MRl. In this respect the current scenario should 
exhibit these characteristics. The MRI operates on an orbital timescale, and so the field would grow in a few tens of milUseconds. It 
may also be amplified by the convective motions described above, § 3.2. The saturation value for the field can be naively estimated as 
that at which its energy density is in equipartition with the gas, B^/Stt w pc,, or when the Alfven speed, = B/ ^Anp is comparable 
with the azimuthal velocity, v^. This gives B « 10'* G. It is not clear at all, however, that the field amphtude will reach such high 
levels, because the magnetic Reynolds number is far beyond its critical value (where diffusion balances dynamo-driven growth), and 
amplification can lead to field expulsion from the convective region, thereby destroying the dynamo (Radler et al. 2002). Precious 
little is known about the growth of magnetic fields at such overcritical levels (Reinhardt & Geppert 2005), and a definitive answer 
will require the self-consistent inclusion of full MHD into the problem at hand (but with a level of resolution which may be well 
above present computational capabilities). It is not clear either that the magnetic shearing instability can generate a mean poloidal 
field as strong as B « lO'* G, since to first order it does not amplify the total magnetic flux. The non-Unear evolution of the instability 
depends sensitively on details of magnetic reconnection, and it has indeed been suggested that this can smooth reversals in the field 
on very small scales, pushing the dominant growing mode to much larger scales (Goodman & Xu 1994). 

It is certainly possible, as shown here, that compact binary mergers do form a neutron disk that is hot enough to be optically thick 
to neutrinos, and convective instability is a direct consequence of the hot nuclear equation of state. A neutron disk is Ukely to be 
convective if the accretion luminosity exceeds 10^"- 10^' erg s"'. Note that even if the accretion luminosity is lower, a hot, massive 
disk (such as those forming in collapsars, Woosley 1993) would undergo a brief period of convection as a result of secular cooling 
(notice that convection is driven by secular neutrino cooling, whereas the MRl is powered by a release of shear kinetic energy). If 
the dense matter rotates roughly at the local Keplerian angular velocity, Q. ~ {GM^^i/r')^!'^, then L^ag is approximately independent 

of radius, and the required poloidal field for a given luminosity is B\s ~ 5q(Mbh/A^o)~' • If a period of convection is a necessary 
step in the formation of a strong, large-scale poloidal field, an acceptable model thus requires that the surrounding torus should not 
completely drain into the hole on too short of a timescale. Whether a torus of given mass survives clearly depends on its thickness 
and stratification, which in turn depends on internal viscous dissipation and neutrino cooling. 

A large amount of differential rotation (as may occur in newborn neutron stars or those in X-ray binaries, and is definitely the 
case in toroidal structures supported mainly by centrifugal forces), combined with short periods, may produce substantial magnetic 
field amplification (Kluzniak & Ruderman 1998; Spruit 1999). The energy transferred to the magnetic field is released in episodic 
outbursts when the buoyancy force aUows the field to rise to the surface of the star or disk. The ampUfication of a magnetic field to 
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such strong values would clearly have important consequences on the evolution and time variability of the disk and its energy output. 
It would probably lead to strong flaring and reconnection events accompanied by the release of large amounts of energy, if the growth 
time for field ampUfication, ?b. is shorter than the accretion timescale, ?acc (otherwise the disk would drain of matter before the field 
had the chance to reach large values; in this case, the survival of a massive, rapidly rotating neutron star as the end-point of binary 
NS merger might be preferred over the prompt formation of a BH). An effective helical dynamo of the a-il type should be favored 
by a low effective viscosity, a, because, as stated in §3.2, the overturn time is ?con ^ 20 ms (this applies only if the disk is not fed 
with matter externally for a time longer than ?acc. otherwise convection would also be able to amplify the magnetic field). 

4.3.3. Nucleosynthesis 

Core collapse SNe and compact object mergers are natural astrophysical sites for the production of heavy elements (Lattimer & 
Schramm 1974; Meyer & Brown 1997; Freiburghaus et al. 1999; Rosswog et al. 1999; Lee 2001). In particular, nucleosynthesis in 
neutrino-driven winds is an issue that may be relevant for iron-group elements, as well as for heavier nuclei through the r-process 
(Woosley & Hoffman 1992). Initial investigations into this matter (Qian & Woosley 1996) determined that the entropy in the outflow 
arising from a newborn neutron star was probably too low to give rise to the r-process efficiently. More recently, this problem has 
been addressed again in the specific case of coUapsar or post-merger accretion disks, based on the results of analytical calculations 
of neutrino cooled disks in one dimension (Pruet, Thompson & Hoffman 2004). If the wind consists of a uniform outflow driven 
from the surface of the disk, the entropy is too low, and essentially only iron-group elements are synthesized, in agreement with the 
earlier results. However, their results also indicate that in a bubble-type outflow (against a background steady wind), where episodic 
expulsion of material from the inner regions takes place, material may be ejected from the disk and preserve its low electron fraction, 
thus allowing the r-process to occur. The convective motions reported here, occurring in the optically thick portion of the disk, would 
represent one way such cells could be transported to the disk surface if they can move fast enough to preserve the neutron excess. 
Since the existence of a convection region is dependent upon the densities reached in the inner disk (so that it becomes opaque), the 
synthesized nuclei (iron group vs. heavier, r-process elements) could be a reflection of its absence/presence. 

We have benefited from many useful discussions and correspondence with U. Geppert, N. Itoh, K. Kohri, P. Kumar, A. MacFadyen, 
P. M^szaros, M. Prakash, M. Rees, T. Thompson, S. Woosley and A. Socrates. Financial support for this work was provided in part by 
CONACyT 36632E (WHL, DP) and by NASA flirough a Chandra Postdoctoral FeUowship award PF3-40028 (ERR). Part of fliis work 
was done during visits to the Institute for Advanced Study (WHL) and Instituto de Astronomia, UNAM (ERR), whose hospitaUty is 
gratefully acknowledged. We thank the anonymous referee for helpful conmients and suggestions on the initial manuscript. 
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Fig. 1. — Color coded logarithmic contours of density (g cm"'), temperature (MeV), electron fraction, optical depth, cooling (erg g"' s"') and magnitude of 
velocity (cm s"') for run a2M at ? = 100 ms. The qualitative change in composition is clearly seen in the contours of Ye when the material becomes optically thick 
(note also the exponential suppression of cooling in the contours of q in this regime). See text for details. 
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Fig. 2. — Density and entropy per baryon, s/k, along the equatorial plane, z = 0, for run a2M at ? = 100 ms, as in Figure 1. Note the logarithimc scale in density 
and the sharp transition and change in behavior from the optically thin to optically thick regimes. Reference power laws for the density profile are indicated. In the 
optically thick region, the entropy per baryon is practically constant. 
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Fig. 3. — Color coded electron fraction in the p—T plane for runs a2M (left column) and alM (right column) at f = 100 ms (top) and t = 200 ms (bottom). 
The solid curving line marks the transition in composition from free nucleons (at high temperatures and low densities) to a particles (at low temperatures and high 

densities). The solid straight line marks the degeneracy threshold, given by kT = 7.7p|(^ MeV. 
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Fig. 4. — Velocity field (left column) and color coded velocity contours (right column) for run a2M at ? = 100 ms (top) and t = 200 ms (bottom), showing clearly 
the small scale circulations present in the disk (the units are cm s"'). The solid white contours correspond to optical depths r = 3/2,6. 
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Fig. 6. — Same as Figure 3, but color coded according to optical depth. 
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Fig. 8. — (a) Toomre parameter QtW = WrCj/TrGS forrun a2M at f = 50 ms and f = 200 ms. (b) Heatin 
(equatorial) disk temperature, Tc, also for run a2M at f = 50 ms. See text for discussion. 
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Fig. 12. — Neutrino luminosity L^, computed from equation (22), for runs alM, a2M and a3M. Note the fairly flat curves and rise, for run alM at early times, 
when the disk is optically thick, and the transition to the optically thin regime at a later time. The curves for runs a2M and a3M exhibit a break approximately on 
the cooling timescale, fcool- 
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Table 1 



Initial conditions for the accretion disks. 



Run'= 


a 








alM 


0.0 


3.85 


0.308 


19,772 


alM 


0.1 


3.85 


0.308 


19,772 


a2M 


0.01 


3.85 


0.308 


19,772 


a3M 


0.001 


3.85 


0.308 


19,772 


aim 


0.1 


3.85 


0.062 


19,772 


a2m 


0.01 


3.85 


0.062 


19,772 


a3m 


0.001 


3.85 


0.062 


19,772 



^Values at the start of the two-dimensional calcula- 
tion of the accretion disk evolution. 

''Number of SPH particles at the start of the two- 
dimensional calculation. 

'^In the run label, the number refers to the value of 
a (with "1" standing for inviscid, with a = 0) and the 
letter to high (M) or low (m) disk mass. 



Table 2 



Accretion disk parameters during the dynamical evolution. 



Run 




pcj/(ergcm ^) 


L,(ergs-') 




£^(MeV)'' 


7'.j/2(ms)'' 




alM 





6-10» 


2-10^^ 


1.5-10" 


8 


>140' 


0.308 


alM 


7 


2-I0-" 


2- 10=' 


2-10" 


8 


60 


0.098 


a2M 


0.7 


2-10" 


6-10" 


2-10" 


8 


150 


0.2 


a3M 


0.05 


1.6-10" 


2.-10'^ 


1.55-10" 


8 


>140' 


0.287 


aim 


0.7 


4-10'° 


2-10'^ 


5.4-10" 


8 


25 


0.022 


a2m 


0.05 


3-10'" 


i-io'^ 


5.5-10" 


8 


75 


0.045 


a3m 


0.005 


3-10™ 


1.5-10" 


7-10" 


8 


>160'' 


0.059 



'^Values are given at / = 0.2 s. 

""This is the duration of the neutrino emission, measured as the time needed for the disk to release one half of . 

"This is a lower limit, since for a = 0.001 the disk is still in a quasi-steady state at t = 0.2 s and has not yet drained into 
the black hole. 

■"The typical neutrino energy is given at the radius where r„ = 1. 



